Spin-dependent gradient correction for more accurate atomization energies of 

molecules 
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We discuss, simplify, and improve the spin-dependent correction of L.A. Constantin et al, Phys. 
Rev. B 84, 233103, for atomization energies, and develop a density parameter of the form v oc 
|Vn|/n 10,/9 , found from the statistical ensemble of one-electron densities. The here constructed 
exchange-correlation generalized gradient approximations (GGAs), named zvPBEsol and zvPBEint, 
show a broad applicability, and a good accuracy for many applications, because these corrected 
functionals significantly improve the atomization and binding energies of molecular systems, without 
worsening the behavior of the original functionals (PBEsol and PBEint) for other properties. This 
spin-dependent correction is also applied to meta-GGA dynamical correlation functionals combined 
with exact-exchange; in this case a significant (about 30%) improvement in atomization energies of 
small molecules is found. 

PACS numbers: 71.10.Ca,71.15.Mb,71.45.Gm 



INTRODUCTION 

Kohn-Sham (KS) density functional theory (DFT) is 
at present one of the most powerful computational tools 
in quantum chemistry, solid-state physics, and materials 
science. Therefore, it is continuously subject to intense 
research. In particular, the development of new improved 
exchange-correlation (XC) functionals [H, which are a 
main ingredient of the KS self-consistent one-particle 
equations Q , and that account for all the quantum many- 
body effects beyond the Hatree approximation, has been 
for many years, and still is, the main topic in this field. 

The construction of XC functionals can be pursued 
in several different ways [S], but in any case a promi- 
nent role is played by the availability of reference sys- 
tems used to target the functional development. In this 
context particularly useful are model systems, which al- 
low to control in a simple way specific features of the 
physics of electronic systems, and thus to reduce the em- 
piricism. Among the others we can mention the uniform 
electron gas [![, which is the basic model for solid-state 
physics and the base of all non-empirical semilocal XC 
functionals 0, Q ; the Airy gas @ and jellium surfaces 
, that describe the physics of simple metal surfaces and 
were used to develop generalized gradient approximations 
(GGA) for solids (AM05 Q, PBEsol Q) and interfaces 
(PBEint [13, El); the hydrogen atom, which is a basic 
model for atomic physics and was successfully employed, 
through the analysis of the XC hole, to construct accu- 
rate meta-GGA approximations I12h14|: the semiclassi- 
cal many-electron neutral atom [lot [l6[ that was used 
to obtain a modified second-order gradient expansion for 
both the kinetic 17 1 and exchange energy [18|, which 
are at the base of accurate XC GGA for molecular sys- 
tems (APBE [HHl), as well kinetic energy GGAs con- 



structed for the embedding theory of weakly-interacting 
molecular systems 2l| . 



Recently, we also proposed [22| a constraint for atom- 
ization energies at the GGA level, based on a model sys- 
tem constituted by a statistical ensemble of one-electron 
densities. In particular, we considered the reference 
one-electron hydrogen (H), Gaussian (G), and cuspless 
(C) densities (atomic units are used throughout, i.e., 
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and we showed that the accuracy of a GGA XC functional 
to describe cohesive properties is directly related to the 
information-entropy-like function [231 ] 

I[GGA] = - Y,PiRi MPi R i), i = H,G, and C, (2) 



where R — \AE XC / E% xact \ is the relative absolute error 
of the XC energy of G, H, or C (for a given GGA), and 
Pi = 1/3 are probabilities (weights) of G, H, and C ap- 
pearance in a molecular bond (using the assumption of 
equ-probability) . Note that this constraint for atomiza- 
tion energies, was derived from an empirical observation 
relating errors in the model one-electron densities to er- 
rors in the atomization energies, of popular GGAs (see 
Fig. 2 of Rcf. [22]); and from the physical explanation 
that one-electron densities are simple models for simple 
bonding regions, where iso-orbital regime can be signif- 
icant, (see Fig. 1 of Ref. [22j, and the corresponding 
discussion.) 

This constraint was used to improve the description of 
cohesive energies of finite-size and periodic systems by 
GGA functionals of PBE-like form [4( through the use of 



2 



a spin-dependent XC correction based on the ansatz [2^ 



d GA («) = e^ A F GGA {s)+e^ A +f a {cj>,t)H 
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where e^ DA and e^ DA are the LDA exchange and cor- 
relation energies per particle; F GGA (s) is the exchange 
enhancement factor; H GGA is the PBE-like gradient cor- 
rection for the correlation; <j) = ((1 + C) 2/3 + (1 - C) 2/3 ) /2 
is a spin-scaling factor that lays between 1 (for spin- 
unpolarized systems) and 2 -1 / 3 = 0.7937 (for fully spin 
polarized systems) and £ = (n-f — n i)/ n is the relative 
spin polarization; s and t are the reduced gradients for 
exchange and correlation respectively [24j, |25[ 



|Vn|/2fcF n, t = |Vn|/2fc s 0n , 



(4) 



with fcp = (37r 2 n) 1/ ' 3 being the local Fermi wavevector, 
and k s = (AUf /n) 1 / 2 being the Thomas-Fermi screening 
wavevector; finally the following spin-dependent correc- 
tion has been selected 12211: 



(5) 



which satisfies important physical constraints [22j. The 
parameter a was obtained by the minimization of 
I[GGA], yielding the zPBEsol and zPBEint GGA func- 
tional 0, having GGA=PBEsol, a = 4.8 and 
GGA=PBEint, a = 2.4, respectively. 

We underline that, as for most DFT functionals 
the ansatzes of Eqs. (3) and (5) were not derived 
directly from exact equations/conditions, but are con- 
structed with the aim to fulfill all the exact constraints 
for the original GGA and additionally, the atomization 
energy constraint, derived in Ref. [22|]. Thus, the spin- 
dependent correction can be viewed as an ad hoc correc- 
tion for the whole XC functional, that can be applied to 
GGAs constructed for solids (see Fig. 3 of Ref. [22]). 

By construction the functionals keep unchanged the 
behavior of the original functionals (PBEsol and PBEint) 
for spin-unpolarized systems, and they significantly im- 
prove the atomization energies of molecules and solids 
[22j . In fact, the ansatz of Eqs. (O and ([5]) preserves 
many important constraints of the original GGAs, in- 
cluding the spin-scaling behavior of the XC functional 
in the slowly- varying density limit 26143 1|, and in addi- 
tion enforce the atomization energy constraint. Note also 
that the spin-dependent correction of Eq. ([5]) only acts 
in the rapidly-varying density regime, where the exact 
spin-dependence of the correlation is not known and the 
PBE-like functionals assume the spin dependence of the 
second-order gradient expansion [j, [5(3] • 

Despite their good behavior, however the zPBEint 
and zPBEsol functionals display some non-negligible for- 
mal drawbacks related to the form used for the spin- 
dependent correction (Eq. (0). In particular, because 
the correction is aimed at improving the bonding descrip- 
tion, it should act only in the valence regions of the elec- 
tron density. On the other hand, as it will be discussed 



in this paper, it turns out that in atomic inter-shell re- 
gions t can be large so that, even if is close to I (as 
in the whole core region), the spin-dependent correction 
may be (inappropriately) important. 

In this paper, we consider in more detail the con- 
struction proposed in Ref. [22j and develop a new spin- 
dependent ansatz, beyond the form of Eq. ([5]), which is 
able to implement the constraint of our statistical set of 
one-particle densities removing the previous limitation. 
To this end we introduce a new density parameter more 
convenient for atomization energies, and additional con- 
straints from partially spin-polarized one electron densi- 
ties. Moreover, using the new ansatz we demonstrate the 
utility of a spin-dependent correction beyond the GGA 
level, applying our spin-dependent factor to improve the 
compatibility of meta-GGA correlation functionals with 
the exact-exchange (EXX). 

This article is organized as follow: in Section n we pro- 
pose an improved ansatz for the spin-dependent correc- 
tion, which is based on a new reduced density parameter; 
in Section n we show the results of several GGA function- 
als bearing the spin-dependent correction, for atomiza- 
tion and binding energies, as well for other properties; in 
Section n we develop spin-dependent corrections for meta- 
GGA correlation functionals, that allows a better com- 
patibility with the exact exchange, improving the atom- 
ization energies by about 30%, and finally, in Section „ 
we summarize our conclusions. 



IMPROVED ANSATZ FOR THE 
SPIN-DEPENDENT CORRECTION 

Let us consider, instead of / Q (0, t) a more general form 
for the spin-dependent correction factor to be used in Eq. 
©: 



-oco 3 |C|' 
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where a > and uj > are constants which will be 
fixed from exact constraints, and £ is the relative spin 
polarization (£ — I for fully-spin-polarized systems, and 
C = for spin-unpolarized systems). Here, v is a gener- 
alized spin-independent density parameter for which we 
use the general ansatz 
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with x a parameter to be fixed from one-electron density 
analysis. Note that the form of Eq. ([7]) assures that 
v can describe the density regime with good flexibility. 
In particular, when x — we have v — t(f>, while for 
x = —1/2 we find v oc s. 

The ansatz in Eq. (j6]) satisfies the following exact prop- 
erties (as the ansatz in Eq. (0: 

(i) At a constant density v = 0, thus / Q , W (C,0) = I, 
and Eq. ([3]) correctly recovers the LDA behavior. 



(ii) For a slowly-varying density, v is small and thus 
/q,u(Ci v ) — * 1 + C(|Vn| 3 ), i.e. with a third power in the 
gradient of the density. In this way the corrected GGA 
of Eq. ([3]), performs exactly as the original PBE-likc 
correlation one 

(iii) For any spin-unpolarized system (£ = every- 
where) or regions of space where C = 0, there is no spin- 
dependent correction (/ a , w (0,i;) = 1), and thus the cor- 
rected GGA of Eq. ([3]) recovers the original PBE-like 
correlation. 

Thus XC-energy density in Eq. ([3]) with the spin- 
dependent correction in Eq. (j6]) satisfies the same phys- 
ical constrains as the original PBEsol (PBEint). In ad- 
dition the spin-dependent correction in Eq. ([5]) is con- 
structed so that fa,u(Ci v ) < 1 onr y m rapidly- varying 
spin-polarized density regime (£ ^ and v > 1), where 
no exact constraints are know. In this regime the to- 
tal correlation energy density of our corrected functional 
(i.e. 6c SDA + fH GGA ) will be larger (more negative) than 
the original PBE-like correlation and can thus partially 
compensate the underestimated contribution of PBEsol 
(PBEint) exchange functional 221 ] . 

Concerning the relation between the new ansatz in Eq. 
([6|) and the one in Eq. ([5]), we note that the latter 
can be rewritten as (f) at = e ln ^ at . For small spin- 
polarizations (C — > 0, <f) — > 1), we have ln(0) — > — |C| 2 /9, 
and v — > t(f> — s> t for x = and <j) — > 1. Thus, the new 
ansatz in Eq. ([6]) is equivalent to the one of Eq. ([5]) with 
lu = 2 and x — {a is a constant). This equivalence 
is trivially valid also for fully polarized systems (as the 
model one-electron densities in the statistical ensemble) : 
in fact in this case <j) — 2~ 1 / 3 and Q — ±1 so that ln(</>) is 
just a constant, = 1 and again v is proportional to t 
for x = 0. 

The main advantage of Eq. (J6j) with respect to Eq. §5§ 
is that the dependence of the correction factor / from 
the spin- and spatial-properties of the density is decou- 
pled and controlled respectively by the u> factor and the 
density parameter v. Therefore, a fine tuning of the cor- 
rection is possible to make it relevant only in valence and 
bonding regions. 

To fix the parameters in Eqs. ([6]) and (0, we note 
that for fully-spin-polarized systems f(v, C = 1) = e~ av , 
thus the parameters a and x can be found from the 
one-electron-densities statistical ensemble, by minimiza- 
tion of the information-entropy-like function. To perform 
such minimization, the form of the GGA functional en- 



tering Eq. ([3]) must be fixed. As in Ref. [22j, we consider 
in the present case GGA=PBEint and GGA=PBEsol, 
obtaining x = 1/6 and a = 1.0 and a = 1.8 for PBEint 
and PBEsol, respectively. 

The minimization procedure for x is also indicated in 
Fig. [U where we plot the function 
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FIG. 1: J(x) (defined in Eq. ©) for PBE, PBEsol, PBEint, 
and Wu-Cohen GGA 



for several popular functionals. The functionals con- 
structed for solids and interfaces (PBEsol, PBEint, and 
Wu-Cohen [32|) show a similar trend with x: the value 
of J is large for x — —0.5, that represents the spin cor- 
rection e~~ as ; for x — 0, corresponding to / = e~ at , i.e. 
the previous ansatz used in zPBEsol and zPBEint, it is 
significantly smaller and close to the absolute minimum; 
at x = 1/6, all three curves minimize the function J(x). 
For functionals more accurate for atomization energies, 
as PBE, the spin-dependent correction is not important 
(and thus J(x) is almost constant), as shown in Fig. [TJ 

We note at this point that the ansatz of Eq. ([7]) used 
to define the density parameter v and the value x = 1/6, 
obtained from the minimization procedure, are impor- 
tant outcomes that go beyond the optimization of the 
parameters of the model given by Eq. ([6]). Thus, the 
minimization of J(x) allows us to define the density pa- 
rameter 
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(9) 



J(x) = min/fa, x] , 



(8) 



with k v — 2(3/(47r 4 )) 1 / 18 n 1 / 9 , being indeed a physi- 
cally meaningful quantity to describe density variations 
in valence and bonding regions. In fact, unlike s and t 
which were derived in the slowly- varying density limit, 
v was derived from the statistical constraint involving 
the one-electron densities, which are characterized by a 
rapidly-varying regime and connected to cohesive ener- 
gies of different systems. Moreover, Fig. [2] shows that, 
in agreement with the indications of wave vector anal- 
ysis on jellium surface energies 33, lllj], k v and the ra- 
tio k v /2kp provide small values of the wave vector in 
the high-density regions (small r s ) and relatively large 
values (~ 1) in valence and bonding regions (r s > 3), 
thus granting a good sensitivity for density variations in 
a broad range of r s values, and a superior description (e.g. 
with respect to k s ) of high wave- vector contributions in 
low-density valence and bonding regions. 
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FIG. 3: E C (Q - E C (C = 0) (hartree) versus ( for the one- 
electron Gaussian densities with uniform spin-polarization £ 
(see Eq. JTOj), for PBEsol, zPBEsol, and zvPBEsol. 



FIG. 2: Dependence of wave vectors Uf, ks, and k v on the 
Seitz radius r a (a.u.). The bottom panel reports the ratio of 
wavevectors k s and k v with 2kp. 



Wavevector analyses for the correlation energies have 
been perfor med on ly for very simple s yst ems, as jellium 
surfaces 
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33|-|35j and Hooke's atom [36j. For such sys- 
tems, the random phase approximation (RPA) fails at 
large wavevectors (k/2kp > 0.8), whereas the semilocal 
functional performs good (see Fig. 9 of Ref. |35|), but 
there is still room for improvement 1361 . We also recall 
that the Langreth and Mehl GGA |24|, which was one 
of the first semilocal functionals beyond the gradient ex- 
pansion, is based on the wavevector decomposition of the 
XC energy. k v , which was found from one-electron anal- 
ysis (see Fig. 1), and agrees well with the Thomas-Fermi 
wave vector (k s ) (see Fig. 2), is more sensitive in rapidly- 
varying regions. Thus, the density parameter v and the 
local wave vector k v , measuring the range of density vari- 
ations important for bonding properties, can be regarded 
as important ingredients in functional development and 
further work shall be planned to fully assess their impor- 
tance in this context (beyond the present spin-dependent 
correction). 

The factor of 3 that enters the v expression was cho- 
sen only for convenience, i.e. for any x, the parameter a 
which minimizes I[GGA] to be between and 2. In fact, 
the minimization of J(x) can define the density param- 
eter v up to a constant, that later may be fixed from a 
theory of rapidly-varying density regions. Nevertheless, 
the here defined k v performs similarly with k s in a broad 
range of density regimes (see Fig. 2). 

The parameter u) could not be fixed from a minimiza- 
tion of the statistical-cntropy-likc function, because for 
fully spin-polarized systems the ansatz of Eq. ((6]) does 
not depend on ui. We must therefore introduce additional 
constraints based on partially polarized model systems. 
To this end we consider the uniformly spin-polarized 
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(10) 



that are models for the electronic localization in Wiener 
crystals and in bonding regions of organic molecules [37j , 
and were already employed to fix constraints in the de- 
velopment of meta-GGA functionals [IHEI, 28, 39j| . The 
following conditions are imposed: 
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forC > 0.7, 



(11) 



where the first condition shall guarantee that the spin- 
dependent correction doesn't modifies the original func- 
tional in the core region of spin-polarized atoms (where 
< C < 0.3) while the second one requires the GGAs us- 
ing the new spin-dependent correction to be the same as 
the previous zGGAs (i.e. zPBEsol or zPBEint), whose 
behavior at large £ had been proved to be remarkably 
good [22 1 . 

We find that for both cases considered here 
(GGA=PBEint and PBEsol) the optimum value of w is 
9/2. Thus, we can finally define two non-empirical func- 
tionals making use of the spin-dependent correction of 
Eq. ([6]) with the density parameter of Eq. ([9]), hav- 
ing u = 9/2 and a = 1.0 (zvPBEint) and a = 1.8 
(zvPBEsol), respectively. 

The plot of E C (C) - E c (( = 0) versus C for 
the one-electron Gaussian densities with uniform spin- 
polarization (Eq. (|10[) ) is shown in Fig. [3] for PBEsol, 
zPBEsol, and zvPBEsol. Similar results (not reported) 
are found for the PBEint case. Here we consider the only 
correlation because the exchange functional is the same 
in the considered functional. By construction, zvPBEsol 
satisfies very well the constraints of Eq. (fTT|) . Note 
that zPBEsol, based on the old ansatz shows a pro- 
nounced C-dependence already for small values of the 
spin-polarization. 
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r (bohr) 

FIG. 4: Upper panel: Spin-dependent correction factors 
f(v,0 of zPBEsol (Eq. ©) and zvPBEsol (Eq. ©), versus 
radial distance r in N atom; Lower panel: The deviations from 
PBEsol of zPBEsol and zvPBEsol, 4nr 2 n(4 unc - e PBEso1 ) 
(a.u.) versus r in N atom. In the inset, is shown 4-irr 2 ne PBEs ° l 
(a.u.) versus r. 



In Fig. |4]we also plot the behavior of the PBEsol spin- 
dependent correction, and the deviation from PBEsol of 
the correlation energy density, versus the radial distance 
in N atom. As expected, the zvPBEsol improves over 
zPBEsol in the core region of the atom, being smooth 
and very close to the original PBEsol curve, thanks to 
its smaller derivative around £ = 0. In the valence and 
tail regions, on the other hand, zPBEsol and zvPBEsol 
perform rather similarly (by construction) and decrease 
significantly the correlation with respect to PBEsol. 

Finally, we stress that the corrected functionals 
(zvPBEsol and zvPBEint) fulfill all the exact constraints 
of the original functionals and additionally satisfy the 
new XC constraint for atomization energies. The Lieb- 
Oxford bound [5HI1 E x > E xc > 2.27E% DA , is lo- 
cally recovered by the exchange part of the function- 
als, and is globally recovered by the full XC function- 
als (in the same manner as PBE and PBEsol), because 
Kc Es ° l Z E z Z PBEso1 > E P C BE (see Table S12 in the sup- 
plementary material). 



RESULTS OF ZVPBESOL AND ZVPBEINT 
CALCULATIONS 

In this section we report the results of the applica- 
tion of the zvPBEint and zvPBEsol functionals to the 
calculation of the atomization energies and other prop- 
erties of several systems. In detail, we consider the at- 
omization energies of organic molecules (AE6 [HI and 
W4 1441 tests), transition-metal complexes (TM10AE 
test SI), and small gold clusters (AUnAE test [46|). 



Note that for all these systems, especially for those of 
the TM10AE test, the atomization/interaction energies 
are obtained as energy differences between species with 
different spin-polarization. Thus, they constitute a nat- 
ural test for the functionals developed in the previous 
section which carry a spin-dependent correction aiming 
at improving cohesive and interaction energies. 

In addition to the above mentioned properties, we con- 
sidered also different properties of spin-polarized systems 
in order to verify, in analogy with Ref. 0, that the 
spin-dependent correction does not worsen significantly 
properties other than cohesive ones, for which it was 
designed. Namely, we considered the binding energies 
of small organic molecules to gold microclusters (small- 
interface (SI7) test 48]), as a model for spin-polarized 
hybrid interfaces, equilibrium bond lengths of open-shell 
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organic molecules (BL9 test |22j) and transition- metal 
complexes (TM10BL test 1220 . kinetic properties of small 
organic reactions (K9 test |47|), the spin-states calcula- 
tions of Table II of Ref. |22|, collected in the AE$ test 
(iron [4!| and cobalt complexes), and the mean ab- 
solute errors on the total XC energies of several spin- 
polarized atoms and molecules 51] (E5* oms and E™ c o1 ). 
The results are reported in Table Q] as well as in Ref. 
All calculations were performed with a development 
version of the TURBOMOLE program package [52j us- 
ing the def2-TZVPP basis-set [53[ and full-self-consistent 
densities. 

Inspection of TableUshows that all functionals bearing 
the spin-dependent correction (zPBEsol, zvPBEsol, zP- 
BEint, and zvPBEint) systematically and significantly 
improve over the original functionals, for all tests con- 
cerning atomization and interaction energies. On the 
other hand, z- and zv-type functionals perform rather 
similarly, although generally a small improvement is ob- 
tained with the zv-construction with respect to the cor- 
responding z-type functionals (except for the W4 test). 
We highlight also that zvPBEsol and zvPBEint have a 
similar performance as PBE for organic molecules atom- 
ization energies, while zvPBEint is among the best GGAs 
for the TM10AE test, the gold clusters (see Fig. 2 of 
Ref. 46), and the hybrid interfaces [48l |. Finally, we re- 
call that both the z- and the zv-functionals preserve for 
the spin-unpolarized cases the performance of the origi- 
nal PBEsol and PBEint functionals, so that they are very 
accurate for properties like lattice constants of paramag- 
netic solids [2fJ , geometries of closed-shell molecules [20( 
and description of large metal clusters [46| . 

Concerning the non-cohesive properties of open-shell 
systems (equilibrium bond lengths, kinetic properties, 
energy differences between spin states, and absolute XC 
energies of atoms and molecules) we see that the z- and 
zv-functionals have similar accuracy and are often better 
than the original functional. Concerning the comparison 
between the z- and the zv-functionals, the former are of- 
ten better, but the differences with respect zv-functional 
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TABLE I: Mean absolute errors (MAEs) for different tests. (The results of the AEs are expressed as mean absolute relative 
errors). For each group of functionals (type-sol and type-int) the best result is highlighted in bold. The PBE and APBE results 
are also shown for reference. Full results are reported in Ref. IHll 



Test 


units 


PBEsol 


zPBEsol 


zvPBEsol 


PBEint 


zPBEint 


zvPBEint 


PBE 


APBE 






Atomization and binding energies 












AE6 (organic molecules) 


kcal/mol 


34.90 


15.72 


14.21 


24.78 


14.62 


14.06 


14.50 


7.93 


W4 (organic molecules) 


kcal/mol 


21.45 


12.64 


13.72 


15.57 


11.59 


12.32 


10.75 


8.56 


TM10AE (transition metals) 


kcal/mol 


18.29 


10.70 


9.24 


15.58 


10.63 


8.50 


13.47 


8.69 


AUnAE (gold clusters) 


kcal / mol 


4.37 


1.97 


2.43 


2.22 


1.01 


0.85 


0.31 


1 .78 


SI7 (hybrid interfaces) 


kcal/mol 


3.92 


2.88 


3.28 


2.80 


2.72 


2.50 


3.69 


5.80 








Other 


properties 












BL9 (organic molecules) 


mA 


15 


16 


15 


17 


15 


10 


15 


15 


TM10BL (transition metals) 


niA 


19 


17 


18 


18 


17 


17 


23 


26 


K9 (kinetics) 


kcal/mol 


10.83 


9.90 


10.42 


9.17 


8.73 


9.03 


7.51 


6.58 


AEs (spin states) 


% 


102 


74 


87 


90 


72 


81 


79 


70 


Ej* oms (atoms) 


mH 


427 


414 


417 


394 


385 


387 


98 


24 


E™° (organic molecules) 


mH 


718 


710 


714 


G61 


657 


658 


146 


83 



are almost negligible. Thus we can conclude that our 
spin-dependent correction, aimed at improving atomiza- 
tion energies, does not introduce significant artifacts that 
might alter the description of the spin-polarized systems. 

In particular, for the AE$ test (energy differences 
between spin states) all GGAs perform badly, with a 
slightly better performance of z-type functionals, that 
in fact is due to their (undesirable) behavior in the core 
of atoms. This test is a hard case for any conventional 
DFT approach as it is related with problems with the 
description of a multireference character by the KS ref- 
erence system and static correlation issues 54]. The 
spin-dependent correction for atomization energies can- 
not be expected therefore to provide any relevant im- 
provement for this property. We also recall that for rel- 
ative energies of different spin-states of metal complexes 
the PBE exchange is known to provide a very poor perfor- 
mance, while much better results are achieved by employ- 
ing the PBE correlation functional in conjunction with 
the OPTX [H|] exchange [H, 

Thus, we believe that the non-empirical zvPBEint, 
which is accurate for energetical and structural properties 
of paramagnetic bulk solids, hybrid interfaces, molecules, 
as well for surface energies of semi-infinite jellium and of 
simple metals [Til ], can have a broad applicability. 

SPIN-DEPENDENT CORRECTION FOR 
META-GGA DYNAMICAL CORRELATION 
FUNCTIONALS 

In this last section we consider the possibility to apply 
the ansatz of Eq. ([6]) also beyond the GGA level. In this 
case however it will not be related to the atomization en- 
ergies through the information-entropy-like function of 
Eq. because the meta-GGAs are one-electron-self- 
correlation free and no information can be gained from 



fully polarized one-electron densities. Instead, we will 
show that the ansatz of Eq. ^ can be usefully employed 
to improve the compatibility of a meta-GGA dynamical 
correlation functional with EXX. For simplicity we con- 
sider the Hartree-Fock non-local exchange: in Ref. Hi 
we showed that similar results can be obtained using the 
local KS EXX. We tested the TPPS meta-GGA corre- 
lation functional [38|, and the recently developed TPSS- 
loc meta-GGA correlation functional [39( , which has the 
same form as TPSS, but uses for its construction instead 
of the PBE GGA correlation, the PBEloc GGA func- 
tional that has the correlation parameter 

P(r s ,t) = 0.0375 + 0.08 t 2 {l-e~^) , (12) 

which ensures a stronger localization of the correlation 
energy density, thus granting a better compatibility with 
exact exchange [39j. 

The idea of the present work is to improve the compat- 
ibility of TPSS/TPSSloc with exact exchange by using 
the spin- dependent correction 

E zvMGGA = J dr n e - aV °\Cr e MGGA ^ ^ 

to provide a further localization of the correlation en- 
ergy density for spin-polarized densities. Here the same 
spin-dependent correction factor of Eq. (J6j> can be used, 
as also in this case we need to modify the original cor- 
relation only in rapidly-varying spin-polarized density 
regime (( ^ and v > 1), thus preserving all the 
other important exact conditions of the original func- 
tional. The rationale beyond this choice is the fact 
that from a physical point of view the localization of 
spin-polarized densities corresponds to a reduction of 
the effective range of same-spin correlation contributions 
e c [n-j-,0] and e c [0,nj, while for opposite-spin contribu- 
tions e c [n^,nj_] — e c [n-j-,0] — e c [0,nj an increase of the 



effective range is obtained. Thus, we can use the ansatz 
of Eq. (p~3|) to (i) fine tune the reiative ranges of the 
same- and opposite-spin correiation contributions, that 
in TPSS (and PBE) were determined from LDA and/or 
jcllium constraint s l30l . l38j ] . so probabiy overestimate the 
same-spin range [14fl with respect to the opposite-spin 
one; (ii) cut the longer-range part of same-spin contri- 
butions that, when exact exchange is used, are no more 
needed to compensate for the too short range of semilocal 
exchange. 

To build the spin-dependent correction factor and find 
appropriate values for the parameters a and ui, we use the 
uniformly spin-polarized Gaussian densities of Eq. (|10p 
and impose the following constraints: 

(i) At small relative spin-polarizations we must recover 
the original functional, hence we must require 



dE. 



zvMGGA 







C<0.3 



(14) 



This condition, in analogy with the GGA case assures 
that in core regions, where closed-shell configurations 
dominate and (, is small, the functional mimics the orig- 
inal behavior of the uncorrected functional. Similar 
constrains has been used for the (rev)TPSS functionals 
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HQ. 

For 



(ii) For arbitrary values of the relative spin- 
polarization the XC functional must approach as close 
as possible the ^-dependence of the exact XC functional, 
E XC (C)- The latter is known to be constant for one- 
electron densities and for different values of £ [H, l59j |. 
i.e. 



E xc (() = E xc (l) = E x (l), 



(15) 



where E x (l) is the exchange energy for the one-electron 
Gaussian density (the correlation energy vanishes for one- 
electron systems). 

Separating the spin-independent part of the exchange 
and correlation energies (E x and E c respectively) from 
their spin-scaling factors we can rewrite Eq. t|l 5|) as: 



E x (0)f(Q + E c (0)g(O = E x (0)f(l) 



(16) 



As the exact spin-scaling factor for the exchange (for one 
electron systems) is: 



/(C) 



(i + 2 + (i-C) 2 



(17) 



the spin-scaling factor for the correlation (for one electron 
system) must be: 



7(0 = 



2- 



(i + C) 2 + (l - C) 2 



= i-C 2 



(18) 



Therefore, we impose the constraint 
E%° MGGA (Q E z c vMGGA (Q)g(Q, at C > 0.7. (19) 



0.000 



-0.005 



k£ -o.oio 



TPSS 

TPSSloc 

• — • zvTPSS 
•- -• zvTPSSloc 



, TPSS 



-0.015 



-0.020 




-v g(QE c 



FIG. 5: E C (Q (hartree) versus £ of the one-electron Gaussian 
density with uniform spin-polarization £ (see Eq. (JTOj ) , for 
TPSS, TPSSloc, zvTPSS, zvTPSSloc, and ideal E c (z) of Eq. 




TPSS 

TPSSloc 

zvTPSSloc 



r (bohr) 



FIG. 6: e c (a.u.) versus the radial distance r for the N 
atom as obtained from the TPSS, TPSSloc, and zvTPSSloc 
functionals. 



In this way in fact the ^-dependence of the overall 
EXX+C functional (approximately) preserves the exact 
spin-behavior, at Q > 0.7, while performs as the original 
functional at ( < 0.3. We note also that Eq. (fl9|) implies 
that E z c vMGGA (l) = E^ GGA (l) = 0. 

By imposing conditions (i) and (ii) we find to = 9/2 
and a = 8 for meta-GGA=TPSSloc, a = 6 for meta- 
GGA=TPSS. The resulting correlation functionals, to 
be used with full EXX, can be named zvTPSSloc and 
zvTPSS, respectively. In Fig. [5] we show how well the 
zvTPSSloc and zvTPSS dynamical correlation function- 
als satisfy constraints (i) and (ii) in contrast to the orig- 
inal TPSS and TPSSloc functionals, which satisfy con- 
strain (i), but whose ^-dependence differs significantly 
from the desired one for £ > 0.4. 

In Fig. [6] we display instead a plot of the TPPS, 
TPSSloc, and zvTPSSloc correlation energy densities 
of the N atom. The figure shows that, as required, 
the zv-correction does not have any effect inside the 
atomic core and only produces a further localization for 
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the valence (spin-polarized) density, showing the corre- 
spondence of the imposed constraints for one-electron 
uniformly-polarized Gaussian densities with the physical 
requirements for the functional. 

At this point, we need to explain the differences be- 
tween the spin-dependent corrections at XC GGA level, 
and at EXX + meta-GGA level. The exchange function- 
als constructed for solids and hybrid interfaces (xPBEsol 
and xPBEint), have localized holes [35[ and thus they 
are not accurate for atoms, overestimating the atomic 
exchange energies. So, the spin-dependent correction, 
based on the H, G, and C ensemble, improves the XC 
energies of one-electron systems, by making a neces sary 
derealization of the spin-dependent correlation hole [22| . 
On the other hand, at the EXX+mcta-GGA level, one 
electron densities are exactly described, and because the 
EXX hole is delocalized, the dynamical correlation part 
should be more localized. 

To test the zvTPSS and zvTPSSloc functionals we ap- 
plied them to compute atomization energies (AE6 test 
[43(|), barrier heights (BH6 test [47J ) and reaction kinetic 
properties (K9 test[47() of small organic systems. More 
extensive tests on larger systems are instead not con- 
sidered, because in our treatment we do not account for 
non-dynamical correlation, so that only energy properties 
of simple small systems at equilibrium geometry can be 
computed at present. We recall in addition that for spin- 
unpolarized systems the zvTPSS and zvTPSSloc perform 
exactly as the original TPSS and TPSSloc functionals, so 
we did not consider tests involving only spin-unpolarized 
species (e.g. non-bonded interactions). 

The results are reported in Table [TTJ From an inspec- 
tion of the table appears clearly that the zv-corrected cor- 
relation functionals have a significantly higher compati- 
bility with exact exchange and provide much improved 
results with respect to the original TPSS and TPSS- 
loc functionals. In particular, the zvTPSSloc functional, 
joining the spatial localization of the correlation energy 
density and the spin-dependent correction, provides in 
all cases the smallest errors, with MAEs of atomization 
energies comparable with those of PBE and very small 
errors for barrier heights and the K9 test. 



TABLE II: Mean absolute errors (kcal/mol) on atomization 
energies and kinetics of small organic systems, as resulting 
from calculations using Hartree-Fock exchange and different 
meta-GGA correlation functionals. The best value for each 
line is shown in boldface. Full results are reported in Ref. l5ll . 



Test set 


TPSS 


TPSSloc 


zvTPSS 


zvTPSSloc 


AE6 


29.1 


25.5 


20.6 


17.3 


BH6 


4.7 


3.9 


4.6 


3.4 


K9 


5.6 


4.3 


4.3 


3.6 



CONCLUSIONS 

In conclusion, we have discussed in detail the GGA 
spin-dependent correction of Ref. |22|, for atomization 
energies. From an analysis the one-electron statistical 
ensemble (H, G, and C), we have found a density pa- 
rameter v oc |Vn|/n 10 / 9 for the valence and tail density 
regions, whose wave vector k v cx n 1 / 9 , seems the next 
term is a wavevector series, where the first two terms are 
kp oc n 1 / 3 (compatible with exchange), and k s cx n 1 / 6 
(compatible with Yukawa interaction 4] in the slowly- 
varying high-density limit). 

By using the density parameter v, we propose a sim- 
pler GGA spin-dependent correction, constructed to be 
relevant only in the valence and tail regions, where the 
spin polarization is important. The here proposed GGAs 
(zvPBEsol and zvPBEint) systematically improve the at- 
omization and binding energies of molecular systems with 
respect PBEsol and PBEint, preserving the accuracy of 
the original functionals for other properties, and thus 
achieving a broader applicability. 

The results of Table I, show that the zvPBEint func- 
tional performs remarkably well for energetical and struc- 
tural properties of transition metals, being more accu- 
rate than the APBE GGA reference, and thus they can 
be considered in applications of transition metal chem- 
istry, where popular functionals (including hybrids) give 
a modest behavior. 

The z- and zv- functionals, perform similarly for atom- 
ization and binding energies, but better results are found 
from zv- type, especially for the PBEint case. On the 
other hand, for other properties, the zv- functionals are 
closer to the original ones (by construction). The con- 
struction of non-empirical GGAs with broad applicabil- 
ity (i.e., accurate for molecules, bulk solids and surfaces), 
is of great theoretical and practical interest. Recent work 
proved that non-empirical GGAs can not be both accu- 
rate for atoms and solids despite highly empirical 
GGA functionals demonstrated a good accuracy for a 
broad range of problems 6l|. Thus, the development of 
non-empirical GGA functionals with broad applicability, 
is a theoretical challenge in DFT. In this work we have 
shown that the spin-dependent correction (applied to sat- 
isfy the statistical condition for atomization energies, and 
other physical conditions derived from model systems), 
may be one path through solving this challenge. 

We have also applied the spin-dependent correction to 
meta-GGAs dynamical correlation functionals, showing a 
significantly better compatibility with the EXX. In fact, 
the zv-corrected meta-GGAs improve more than 30% 
over the original ones, for atomization energies. Thus, we 
believe that the zvTPSSloc meta-GGA correlation func- 
tional, that is accurate for jellium surfaces, and Hooke's 
atom at any frequency [3^ | , can be a good starting point 
in the development of more accurate (and non-empirical) 
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hyper-GGAs [62J, 16 



Finally, in this paper, we have shown that the spin- 
dependent correction of Ref. l22l which was simplified 
and improved in this work, is a powerful tool at GGA 
and meta-GGA levels. Thus, we believe that it can be 
used also to functionals that describe the dimensional 
crossover (from 3D to 2D) of the XC energy [Hj, in or- 
der to improve the atomization and binding energies of 
molecular systems under a 2D confinement. 
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